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ABSTRACT 

We have apphed the torus fitting procedure described in Ng & Romani (2004) 
to PWNe observations in the Chandra data archive. This study provides quanti- 
tative measurement of the PWN geometry and we characterize the uncertainties 
in the fits, with statistical errors coming from the fit uncertainties and systematic 
errors estimated by varying the assumed fitting model. The symmetry axis \E' of 
the PWN are generally well determined, and highly model- independent. We often 
derive a robust value for the spin inclination C,. We briefiy discuss the utility of 
these results in comparison with new radio and high energy pulse measurements. 

Subject headings: neutron- stars: rotation-stars: winds, outfiows 

1. Introduction 

One of the greatest success of the Chandra X-ray Observatory (CXO) is the discovery 
of equatorial tori and polar jet structures in many pulsar wind nebula (PWN) systems. It is 



now b elieve that these features are common among young neutron stars. In the lRees fc Gunn 



( I1974J ) picture, when the highly relativistic pulsar wind decelerates in the external medium, 

aci 

(. 



a termination shock is formed at a characteristic scale 

1/2 



ext 



where E is the pulsar spin-down power and rj is the filling factor. In general, if the pulsar 
is subsonic in the ambient medium, i.e. Pcxt > -Pram = 6 x lO^^^nvj g cm^^s^^ for a pulsar 
speed lO'^fy cm s~^ through a density of nnip cm"^, a toroidal shock structure is expected; 
faster objects produce bow shock nebulae. Indeed, many young pulsars still inside their high 
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pressure supernova remnant birth sites do show such toroidal symmetry. The best-known 



exara ple is the PWN around the Crab pulsar as ob served by the CXO (IWeiss. 



copf et al. 



20011). Recently, severa l relativistic MHD models, e.g. iKomissarov fc LyubarskyI (120031 ) and 
Del Zanna et al.l (120061 ) . have shown how such toroidal structure can form if the pulsar wind 
has a latitudinal variation. 



Ng fc Romanil (120041 ) (hereafter iPaper j ) developed a fitting procedure to measure the 



3D orientation of the pulsar wind torus and applied to a few X-ray observations. While this 
simple geometrical model does not capture the fine details of the MHD simulations, is does 
allow one to extract the to rus (and hence pulsar spin) orientation from relatively low signal- 
to- noise ratio (S/N) data. iPaper 11 also gave quantitative estimates for the statistical errors 
arising from Poisson statistics. However the systematic errors due to unmodeled components 
such as jets or background were neglected. For bright objects e.g. the Crab and Vela pulsars, 
the S/N is high and such systematic errors dominate. In this study, we attempt to quantify 
these systematic errors and apply the fitting to more CXO PWN observations, thus providing 
a more comprehensive study. 



2. Pulsar Wind Torus Fitting 
2.1. Brief Review 



In iPaper j . we developed a generic torus fitting procedure to capture the characteristic 
structure of PWNe. Here we give a brief summary of the model: consider a torus in 3D with 
radius r and Gaussian cross section of thickness 6. The emissivity is described by 



Iq cx exp 



-J^ + i^/x'^ + y'^ -r) 
252 



where z' is the direction along the torus axis and the observer line of slight lies in the y'z'- 
plane. As the post-shock flow in the PWNe is expected to be mildly relativistic, the torus 
intensity in the observer's frame is Doppler boosted by 

J oc (1 - n ■ /3)-(^+^)Jo . 

Here we assume a constant flow speed j3=v/c radially from the pulsar in the equatorial plane. 
Of course this does not captu re the complex latitud inal structure and speed variations seen 
in the MHD simulations (e.g. iDel Zanna et al.ll2006l ). but it does give a characteristic speed 
and captures the gross Doppler boosting. During the fitting we leave /3 as a free parameter, 
but assume a fixed photon spectral index P = 1.5 in the rest frame, which is a typical value 
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for PWNe. We have also tried other va lues of F, e.g. F = 1.2 and F = 1.9 as report ed for 
the Vela and Crab pulsars respectively (jKargaltsev &: Pavlov 2004 : Mori et al. 2004). The 
geometrical parameters are quite insensitive to F. The phenomenological speed /3 is slightly 
affected (with 6(3 ^ — 0.15 5F), so large spectral index variations could induce trends in 
the fit speed. Altogether this gives 



N 



y' sin C 



(27r(5)V 




-(i+r) 



exp 



252 



where is the total number of counts in the torus and ( is the inclination of the torus. 

In order to compare with data, the 3-D torus is projected onto the sky plane. We set up 
the coordinate system for the image plane with y along celestial north, z along the observer 
line of sight and the CCD in the xy-plane. Then the transformation between the frames is 
given by 

a;' \ / - cos \E' - sin 

?/' = I sin \E' cos C — cos \E' cos C sin ( 
z' J \ — sin \E' sin ( cos \I' sin ( cos C 

where is the position angle (N through E) of the projected torus axis. The counts in each 
CCD pixel (x,y) result from the integration of emissivity through the line of sight 



c{x,y) = Y,Hx',y',^') 



The basic torus model is characterized by the parameters \E', (, r, S and f3. In addition, 
a point source and constant background are also included in the fit. The counts in each 
component are determined by the fit, while the total counts in the model is constrained to 
match the image counts. In a few cases where twin tori are observed, there is an extra fitting 
parameter d for the separation. For simplicity we assume a symmetrical offset along the 
torus axis, with the pulsar in the middle. A few other cases have inner and outer tori, at 
different r. 

To estimate the goodness of fit for particular model parameters, a likelihood function is 
defined using Poisson statistics. The probability of observing d^y counts out of the expected 
Cxy counts in pixel {x, y) is given by 



~ d 



xy 



In the high-count limit, we switch to Gaussian statistics for C^y > 20, hence the likelihood 
function passes to the usual statistics. The best-fit parameters are then determined by 



-4- 




Fig. 1. — Illustration of a twin torus, labeled with the fitting parameters. 

the maximum likelihood method, i.e. minimization of the Figure-of-merit function, which is 
defined as the negative log of the likelihood function summed over all image pixels. 

Once the best-fit parameters are obtained, statistical errors in the fit are estimated 
by Monte Carlo simulations. Poisson realizations of the best-fit model are simulated and 
fitted, and the confidence interval of each parameter is obtained from their distributions. To 
estimate the systematic errors, we blank out different regions in the data (i.e. omit those 
pixels from the fit), then re-fit the model by varying only one parameter at a time while 
fixing the others at their best-fit values. The changes in the best-fit parameter values thus 
provide estimates of the uncertainties contributed by un-modeled features in the data. Note 
that in addition to the single parameter systemetic errors computed here, covariance between 
the parameters could also lead to larger errors (e.g. the covariance of f3 with F). 



2.2. Implementation 

In this study we went through the CXO data archive and identified a list of observations 
with apparent toroidal termination shocks in PWNe (Tabled]). There are a few other objects 
with jets and/or possible tori (e.g. PSR B1509— 58), but for these cases closer examinations 
suggest that the pulsar motions maybe trans- or super-sonic with respect to the surrounding 
medium, hence static torus fits are not appropriate. The data are cleaned up with standard 
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reduction procedures using CIAO and the background light curves are examined to filter out 
periods that suffered from strong particle flares. The n we remove the pixel randomization 
for the ACIS observations and apply an algorithm by lMori et al.l (120011 ) . which corrects the 
positions of the split pixel events, to improve the spatial resolution. Finally we select the 
appropriate energy range and binning to optimize the S/N of the PWNe, separating out the 
soft background from the SNR and field stars. To model the point source, a high S/N PSF 
model is simulated with the Chandra Ray Tracer (ChaRT) at the detector source position, 
using the source spectrum reported in literature. 



In the fitting procedure, a simulated annealing package (jPress et al.lll992l ) is employed 
for multidimensional minimization of the Figure-of-merit function, with the initial parame- 
ters chosen by eye. After the best-fit model is obtained, we simulate 500 Monte Carlo Poisson 
realizations and estimate the statistical errors from the parameter distributions accordingly. 
The la confidence interval, which corresponds to 68% confidence level, are reported. To 
estimate the systematic errors, we begin by blanking out the jet regions in the observa- 
tions, since they are the most obvious un-modeled features and are present in many systems. 
As an example, the image of PSR B0540— 69 wit h iet regions remo ved is shown in Fig- 
ure [H The golden section search in one dimension (jPress et al.l 119921 ) is then employed for 
re-fitting the model. This procedure gives systematic error estimates for the Crab, Vela, 
PSRs J2229+6113, B1706-44, J2021+3651 and B0540-69. However for PSRs J1930+1852 
and J0205+6449, the jets are relatively faint and well-separated from the PWNe. There- 
fore, removing the jets alone does not change the parameter values significantly. We believe 
this underestimates the systematic errors and contributions from other un-modeled features 
should also be considered. Therefore, we obtained another estimate by removing the point 
source in the data. The results reported for these two objects are the combination of the 
two systematic error estimates, added in quadrature. 



Table 1. Archival dataset used in this study. 



Pulsar 


SNR 


ObsID 


Instrument 


Exposure (ks) 


J0205+6449 


3C 58 


3832, 4382, 4383 


ACIS 


350 


J0537-6910 


N157B 


2783 


ACIS 


50 


B0540-69 


SNR 0540-69.3 


132, 1735, 1736, 1727, 1738, 1741, 1745 


HRC 


75 


J1124-5916 


G292.0+1.8 


1953 


HRC 


50 


B1800-21 




5590 


ACIS 


30 


J1833-1034 


G21. 5-0.9 


1230, 1554,159, 2873,3700 


ACIS 


56 
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2.3. Individual Objects 

In the following we discuss the fits of some individual objects. 

PSR B0540—69: The combined HRC exposure shows polar jet structure in the PWN. 
In addition, a bright blob SW of the point source is also observed. We have tried blanking 
out the blob in the fit, and found that the result is insensitive. 

PSR J0205+ 6449: Similar to the Crab pulsar, a double tori structure is observed in this 
PWN. However, in this case the outer torus is very diffuse and without a clear boundary. 
In order to have a better constraint in the position angle, we fit the two tori simultaneously 
with the same value of 

PSR Jl 124- 591 6: While the PWN image in the initial CXO HRC observation showed 
noticeable ellipticity (iHughes et al. the elongation was not clear in a subsequent ACIS 



exposure. This could possibly due to the large off-axis angle degrading the PSF. Therefore 
in this study we fitted the HRC image only. A scheduled deep ACIS observation could help 
to resolve this structure. 



PSR J0537-6910: A cometary trail is observed in this PWN system fIChen et al.ll2006t l. 
which is modeled in the fit by a rectangular background region with constant brightness. 
Comparison with the best-fit parameters without the background provides an estimate for 
the systematic errors. 



PSR B1800-21: Although the 30 ks ACIS exposure only detects a total of ~200 counts 
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(0.8-6 keV) for the point source and the PWN, the latter clearly shows extended structure 
with a well-defined symmetry axis. Interesti ngly, the direc t ion o f elongation is in general 
orthogonal to the proper motion claimed in iBrisken et al.l (120061 ). Hence, employing the 
torus fitting to obtain a quantitative measurement of the accuracy is valuable. 

PSR J1833— 1034: As it's host SNR is one of the calibration objects for CXO, the pulsar 
has been observed repeatedly. We combined all the on-axis and in-focus ACIS observations 
in the data archive and obtain a total exposure of 56 ks. 

PSR J0538~2817: In a short 20 ks ACIS exposure, diffuse emission around t he pulsar 
was fi t to a simple equatorial fiow. However, while a deeper 100 ks observation (INg et al. 
20071 ) showed extended emission along a similar axis, it was not unambiguously possible to 
decide whether this was an edge-on equatorial torus or a polar jet. To be conservative, we 
do not include this object in the present analysis. 



3. Results and Discussions 

Table [2] lists the best-fit results with corresponding statistical and systematic errors. 
We report the systematic errors for \E', (, r and the separation only, since they are the most 
interesting parameters. Systematic errors in other parameters have similar fractional values 
when compared to the best-fit parameters. Note that for the faint objects, including PSRs 
J1124— 5916, B1800— 21 and J1833— 1034, the detailed structures are not resolved and the 
large statistical errors dominate over the systematic errors. Therefore only the former are 
reported here. The best- fit mode ls for the new objects in this study are shown in Figure [31 



the rest may be found in iPaper 11 . Figure 2 and the references listed in Table [21 



As shown in the table, the fits for Crab and Vela are dominated by systematic errors 
and the statistical errors are negligible. For the other objects, the two are comparable in 
and r, but the systematic errors for ( are in general much larger. This suggests the fit is 
most robust in obtaining and r, as they are largely determined by the shape of the torus. 
On the other hand, ( depends sensitively on the relative brightness of both sides of the torus; 
thus, removing the jet structure can change the results significantly. 



3.1. Pulsar Jets 

As mentioned in the previous section, polar jets are observed in many of the pulsars in 
this study. In general, the jets' position angle is in reasonable agreement with \E' axis deter- 
mined from the torus geometry. Since these are also relativistic, relative Doppler boosting 



Table 2. Best-fit 'torus' parameters with la uncertainties. The statistical and systematic errors are hsted as the first 
and second error terms respectively. We report the systematic errors for selected objects in (, r and f3 only. 



object 


V( ) 


C(, ) 


r( ) 


) 


p 


Sep { ) 


r^o/ torus i^utsj 


xiei. 


Crab (inner) 


124 ±0.1 ±0.1 


61.3 ±0.1 ±1.1 


15.6 ±0.03 ±0.1 


1.5* 


0-49to:oo6 ±0.006 




■ ■ ■ /l.O X 10^ 


1 


Crab (outer) 


126.31 ±0.03 ±0.11 


63.03lo o3 ± 1.3 


41.33to;°3 ±0-2 


5.9* 


0.55 ± 0.001 ± 0.0004 




■ ■ ■ /l.l X lO'^ 


1 


Vela 


130.63lo °7 ± 0-05 


63.6t°;°^±0.6 


21.25to;o2 ±0.5 


3.0* 


0.4410:00^ ±0.008 


11.61 ±0.03 ±0.4 


■ ■ ■ /1. 3 X 10** 


1 


J1930+1852 


9ll* ± 1.1 


147 ± 3 ± 3 


4.6 ±0.1 ±0.14 


1.1±0.1 


0.62t°;°^±o.l 




1701/602 


1 


J2229+6114 


103±2±1.6 


46 ± 2 ± 6 


9.3 ±0.2 ±0.3 


2.5* 


0.49 ±0.02 ±0.09 




2221/1113 


1 


B1706-44 


163.6 ± 0.7 ± 1.6 


53.3t{ l ±2.9 


3.3l°:o6 ±0.1 ±0.01 


1.0* 


0.70 ±0.01 




2897/1221 


2 


J2021+3651 


45 ±1.3 ±0.6 


79 ± 1 ± 2 


8.0 ±0.2 ±0.4 


1.2* 


0.64 ±0.02 ±0.02 


3-7^0.1 ± 0.07 


235/751 


3 


J0205+6449 (inner) 


90.3t ±0.2 ±0.4 


91.6 ±0.2 ±2.5 


2.08^°;°^ ± 0.01 


0.74* 


0.66lo;gi ±0.009 


5080/39730 


4 


J0205+6449 (outer) 


90.3t ±0.2 ±0.4 


90.56lo;og ±0.01 


16.3lp f ±0.1 


2.0* 


0.88l°;oot ± 0.0007 




5080/29503 


4 


J0537-6910 


131±2±0.9 


92.8^0,1 ±0.5 


4.0l°;i ±0.3 


0.47 ±0.04 


0.86 ±0.02 ±0.08 




102/179 


4 


B0540-69 


144.1 ±0.2 ±0.8 


92.9 ±0.1 ±0.6 


2.35 ± 0.02 ± 0.06 


0.66* 


0.64 ±0.01 ±0.02 




13442/50215 


4 


J1124-5916 


16 ±3 


105 ± 7 


0.9 ±0.1 


'-'■'^-0.05 


15+0-06 
'^■^''-0.04 




390/552 


4 


B1800-21 


44 ±4 


90 ±2 


o 1 +0.4 
■^■^-0.5 


0.74* 






18/88 


4 


J1833-1034 


45 ± 1 




5.7 ± 0.2 


~ 1.2 


0.86 ± 0.01 




304/6541 


00 4 



*- held fixed in the fit. 
t— fitted simultaneously. 

References. — Best-fit parameters and statistical errors are from: (1) iNg fc Romaiiil ll2004l) . (2) iRomani et all l l2005h . (3) iHessels et al] ||2004| '). (4) - this work. 
All systematic errors are new in this work. 
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Fig. 3. — New fitting results in this study. Chandra data (left) compared with the best-fit 
models (right) for PSRs B0540-69, J0205+6449, J1124-5916, J1833-1034, J0537-6910 
and B1800— 21(/e/f to right, top to bottom). 

could in principle provide a second estimate of spin inclination ( and flow speed (3. However, 
two effects complic ate this. First, /3 migh t be quite different in the jets as suggested in the 
MHD simulations (jPel Zanna et al.l 120061 ). Second, the jets are in several cases distinctly 
curved, and repeated exposures of Crab and Vela indicate that this bending in the plane of 
the sky changes with time. Since we expect similar curvature along the line-of-sight, and 
even quite small changes in inclination strongly affect the brightness ratios, interpretation 
of the jet boosting is more difficult than that of the tori. 



3.2. Pulsar Kick Connection 



One direct application of the torus fitting is measurement of the neutron star's spin 
axis orientation. The alignment angle between the pulsar spin and velocity vectors provides 
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insights into the supernova core collapse physics. For a simple picture of a kick acting at 
a fixed spot on the neutron star s urface, a longe r kick duration gives a better alignment 
as a result of rotational averaging. iNg fc Romanil (120071 ) performed simulations of neutrino 
kick models and extracted constraints on the physics parameters including the kick direction, 
magnitude and timescale by comparing torus fitting angles with proper motion estimates and 
other pulsar data. As the sample of X-ray PWNe with modeled torus structure increases, 
the statistical constraints on such kick models should rapidly improve. 



3.3. Radio Polarization 



There are a few objects in ou r list which also have independent estimates for the spin 
geometry from radio polarimetry. iJohnston et al.l (120051 ). in particular, have recently made 
considerable progress in measuring by determining the polarization at the phase of closest 
approach to the magnetic pole and correcting for Faraday rotation to obtain the absolute 
position angle on the sky. Unfortunately, radio pulsar emission is complex and can switch 
between orthogonal modes; this implies a 7r/2 ambiguity in the spin axis orientation. Never- 
theless, careful analysis can give usef ul position angl e s. At present only three of our pulsars 
have a from radio measurements. iJohnston et al.l (120051 ) observed the radio polarization 
for the Vela pul sar and obta i n a P A. of 126?8 ± 1°. For PSR B1 706-44, the res u lt is c on- 
troversial: while IWang et al.l J2006h report a P.A. of 162?0 ± 10?0, Ijohnston et~all J2005h do 
not give a value, since in new observations the phase of closest approach is difficult to locate. 
Finally, iKothes. Reich fc Uyanikerl (120061 ) measure the radio polarization of the PWN in 
PSR J2229-I-6114 and suggest a symmetry axis at P.A.~ 93°. While these three values are 
in reasonable agreement with the \E' from the X-ray symmetry axis, more data are needed 
for a serious comparison. 

Unfortunately, interpretation of the polarimetry appears to be particularly challenging 
for the young pulsars showing X-ray tori (this may be related to their powerful accelera- 
tors in relatively small magnetospheres) . For example, the main radio pulses in the Crab 
appear at the same phases as the high energy emission and may indicate high-altitude emis- 
sion. A classic polarization sweep may only be found for the 'precursor', a low frequency 
component. However, measurement of radio polarization for more objects would be very 
useful. Comparison with the X-ray P.A.s can reveal orthogonal emission in the radio, e.g. 
Helfand. Gotthelf fc Halperru (120011 ). and may help to calibrate the radio polarization tech- 
niques, which may then be applied to older pulsars. PSRs J2020+3651 and B1800— 21 seem 
amenable to such analysis. 



In principle, for a simple dipole geometry it is also possible to measure the magnetic 
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pole impact parameter (3 (from the polarization PA sweep rate at closest approach) and the 
magnetic inclination a (from the overall shape of the polarization position angle sweep), 
when the polarization is well measured over a large range of pulsar phase. In practice, 
while (3 can be fairly well determined, the a estimates are often poor and subject to large 
systematic errors. The Rotating Vector Model (RVM) for such fits was develo ped for Vela, so 



not su rprisingly there are angles for this pulsar. Table 3 gives the values from I Johnston et al. 



( I2OO5I ). However, the sum a + (3 = ( = 43° — 6?5 ~ 37° is in poor agreement with the X-ray 
inclination of 64°. It is difficult to interpret the CXO arc structure seen for Vela with such a 
small inclination. Perhaps the best approach is to adopt C from the PWN structure, when 
available, and (3 from radio polarization. For Vela this would imply a ~ 70°. 



3.4. The 7-ray Connection 

Radio (coherent emission) pulse profiles are often quite complex, and may depend on 
the dense plasma and complex field structure at their relatively low emission altitudes. High 
energy (incoherent emission) profiles should be simpler, especially i f they originate at high 



altitud e where only dipole field components are expected to persist. iRomani fc Yadigaroglu 



(119951 ) found that for the generic class of outer magnetosphere models, the 7-ray pulse is 
usually a double peak, with a peak separation A and phase lag 6 from the radio pulse 
primarily determined by (, although the actual values do depend on a and the gap width, as 



well. Other models, employing emission from both hemispheres (e.g. iDyks fc Rudakll2003l ) 
can have more complex light curves with different phase relationships. However, there are 
some generic predictions of these models: 7-ray peaks are strong and wide when ( ~ 7r/2, 
they converge as ( approaches ~ ir/A, and are weak or absent for small angles. 

In Table |3l we list the observed GeV pulse widths (A^) for Crab, Vela and PSR 
B1706— 44 and tentative widths for two other pulsars discovered after the end of the EGRET 
mission. The fi nal column (Ann) lists approxim ate high energy pulse widths (or ranges) from 



the contours in iRomani fc Yadigaroglul (119951 ) . figure 4. In most cases, a range is allowed - 



independent (e.g. radio) values for P would make width predictions. Note that for Crab and 
Vela, this model tends to predict somewhat narrower pulses than those seen. Also the ten- 
tative detection of PSR J2229+6114 has a wider pulse than predicted for its small (. These 
larger angles may be easier to accommodate in two-pole models. At present there is sufficient 
latitude in the model predictions to allow a wide range of pulse shapes. However, new 7-ray 
observations with AGILE and, especially, GLAST should provide many high quality pulse 
profiles. We expect to have X-ray measurements of ( for many of these young pulsars. If 
radio observations can supply (3, the geometry is determined and comparison with the 7-ray 
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data can powerfully refine (or eliminate) the magnetosphere models. 

4. Conclusion 

In conclusion, we have applied the torus fitting technique to more PWN observations 
in the Chandra data archive and characterized the uncertainties in the fits. This study 
provides a better understanding of the systematic errors, giving quantitative estimates of 
the measurement uncertainties. We argue that these robust position angle and inclination 
( values arc particularly useful for comparison with the radio and high energy pulse data. 
If new observations can fill in more measurements from these energy bands in Table 3, we 
should be able to make substantial progress in understanding the emission zones and viewing 
geometries of young pulsars. 

This work was supported by NASA grant NAGS- 13344 and by Chandra grant AR6-7003 
issued by the Chandra X-Ray Center, which is operated by the Smithsonian Astrophysical 
Observatory for and on behalf of the National Aeronautics Space Administration under 
contract NAS8-03060. 

Facilities: CXO (ACIS, HRC) 
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Table 3. Radio/Gamma- Ray Comparisons 



Pulsar 




C 




Or 


/3fl 






Crab 


124 


61 


7 






0.39 


0-0.35 


Vela 


131 


64 


126.8 ±1 


137 (43)=^ 


-6.5 


0.42 


0.35-0.4 


J1930+1852 


91 


147 













J2229+6114 


103 


46 


~ 93t 






(~0.4)b 


<o.i 


B1706-44 


164 


53 


162 ± 10 






0.25 


0-0.3 


J2021+3651 


45 


79 








(-0.5)- 


0-0.45 


J0205+6449 


90 


92 










0-0.55 


J0537-6910 


131 


93 










0-0.55 


B0540-69 


144 


93 










0-0.55 


J1124-5916 


16 


105 










0-0.45 


B1800-21 


44 


90 






~-6d 




~0.5 


J1833-1034 


45 


85 










0-0.55 



tfr om the polarization structure of the radio nebula dKothes. Reich &: Uvanikej 
l2006h . 

Hjohnston et al.l ||2005| ') 

^Thompson et al. I ll2002h 

' iMcLaughlin fc Corded booj ) 

' ISmits et all l|2006ll 
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